home *** CD-ROM | disk | FTP | other *** search
/ The X-Philes (2nd Revision) / The X-Philes Number 1 (1995).iso / xphiles / hp48hor1 / fractals.src < prev    next >
Text File  |  1991-05-29  |  25KB  |  381 lines

  1. %%HP: T(3)A(R)F(.);
  2. DIR                             
  3. @ FRACTALS, by Dan Ciarniello
  4. @
  5. @ These programs will generate the Mandelbrot set.  To speed up calculations, a
  6. @ couple of properties of the Mandelbrot set are used.
  7. @
  8. @ Since the Mandelbrot set is symmetric about the x-axis, if the x-axis is
  9. @ present, the program will generate the image on only one side and then map
  10. @ the image to the other side of the axis.
  11. @
  12. @ The Mandelbrot set is a closed connected set.  This means that if all points
  13. @ on the border are in the Mandelbrot set, then all points interior to the
  14. @ rectangle are also in the Mandelbrot set.  Ditto if all points are not in the
  15. @ Mandelbrot set.  Thus the program calculates all points around the border of
  16. @ a rectangle.  If all points on the border are the same (in or out), then fill
  17. @ the rectangle the appropriate colour and go to the next rectangle.
  18. @ Otherwise, cut the rectangle in half and check the border of the new
  19. @ rectangle.  This is called the Mariani technique.  Note that loop counters
  20. @ run through pixel coordinates.  These pixel coordinates are converted to user
  21. @ coordinates using PX\->C which uses PPAR to determine conversion factors.
  22.  
  23.   CST { MANDEL SAVE               @ Custom menu.
  24. LOAD { DITH                       @ DITH turns dithering on
  25.     \<< 3 SF                      @ NODITH turns dithering off
  26.     \>> } { NODITH                @                     
  27.     \<< 3 CF                      @                     
  28.     \>> } }                       @                     
  29.   LOAD                            @ LOAD restores a previously saved
  30.     \<< 2 MENU                    @ image and its associated PPAR.
  31. "Enter Filename" ""               @
  32. INPUT OBJ\-> OBJ\->               @
  33. DROP PICT STO                     @
  34. 'PPAR' STO 0 MENU                 @
  35. GRAPH                             @
  36.     \>>                           @
  37.   SAVE                            @ SAVE an image along with its PPAR.
  38.     \<< PPAR PICT RCL             @ PPAR is necessary so that the correct
  39. 2 \->LIST                         @ user units are restored when the
  40. "Enter Filename" {                @ image is restored with LOAD.
  41. "" \Ga } INPUT OBJ\->             @
  42. STO                               @
  43.     \>>                           @
  44.   MANDEL                          @ The main program
  45.     \<< ERASE { # 0d              @ Erase PICT and view it as program
  46. # 0d } PVIEW RCLF                 @ runs; save flag status; start
  47. TICKS \-> c1 c2 f t               @ timing; save initial window corner
  48.       \<< 4 CF 5 CF               @ values.
  49.         IF 'nITTR'                @ Check whether nITTR exists.  If
  50. VTYPE -1 ==                       @ not, create it with a default value
  51.         THEN 100                  @ of 100.
  52. 'nITTR' STO                       @
  53.         END                       @
  54.         IF c2 IM c1               @ Correct the aspect ratio of the
  55. IM - c2 RE c1 RE -                @ image.
  56. 64 * 131 / >                      @
  57.         THEN c2 IM                @ If y2-y1 > (64/131)(x2-x1) then set
  58. c1 IM - 131 * 64 /                @ x-value of upper right corner (x2) to
  59. c1 RE + c2 IM                     @ x1+(y2-y1)*131/64
  60.         ELSE c2 RE                @
  61. DUP c1 RE - 64 *                  @ If y2-y1 < (64/131)(x2-x1) then set
  62. 131 / c1 IM +                     @ y2 to y1+(x2-x1)*64/131
  63.         END R\->C                 @
  64. 'c2' STO                          @
  65.         IF 'PPAR'                 @
  66. VTYPE -1 ==                       @ Check whether PPAR exists.  If not,
  67.         THEN {                    @ create it.                          
  68. (-2,-1) (2.09375,1)               @                                     
  69. X 0 (100,100) TRUTH               @                                     
  70. Y } 'PPAR' STO                    @
  71.         END PPAR 2                @                                     
  72. c2 PUT 'PPAR' STO                 @ Store the window corners in PPAR to
  73. PPAR 1 c1 PUT                     @ set the user coordinates.           
  74. 'PPAR' STO                        @
  75.         IF c1 IM c2               @                                     
  76. IM * 0 <                          @ Determine whether the x-axis is     
  77.         THEN c1 RE                @ present.  If so determine the pixel
  78. 0 R\->C C\->PX c2 RE 0            @ coordinates of the left and right   
  79. R\->C C\->PX \-> lside            @ end.                                
  80. rside                             @                                     
  81.           \<<                     @                                     
  82.             IF c2                 @ Determine whether the x-axis is in  
  83. IM c1 IM ABS >                    @ the upper or lower half of the      
  84.             THEN 0                @ screen.  If its in the lower half,  
  85. 130                               @ then generate the image from top of
  86.               FOR x               @ the screen down to the x-axis.      
  87. x R\->B # 0d 2 \->LIST            @ ITERATE around the first rectangle  
  88. PX\->C ITERATE                    @ starting at the top of the screen.  
  89.               NEXT                @                                     
  90. 1 rside 2 GET B\->R               @                                     
  91.               FOR y               @ ITERATE down right side to x-axis.  
  92. # 130d y R\->B 2                  @                                     
  93. \->LIST PX\->C ITERATE            @
  94.               NEXT                @
  95. 129 0                             @
  96.               FOR x               @ ITERATE along the x-axis from right
  97. x R\->B lside 2 GET 2             @ to left.                            
  98. \->LIST PX\->C ITERATE            @                                     
  99. -1                                @                                     
  100.               STEP                @                                     
  101. lside 2 GET B\->R 1               @                                     
  102.               FOR y               @ ITERATE up the left side of the     
  103. # 0d y R\->B 2 \->LIST            @ screen from the x-axis to the top   
  104. PX\->C ITERATE -1                 @ of the screen.                      
  105.               STEP                @                                     
  106. { # 0d # 0d } rside               @
  107. MARIANI 63 lside 2                @ Pass the upper left and lower right
  108. GET B\->R                         @ pixel coordinates of the rectangle  
  109.               FOR y               @ to MARIANI.                         
  110. # 0d y R\->B 2 \->LIST            @                                     
  111. PICT OVER PX\->C CONJ             @ Map the area above the x-axis to    
  112. C\->PX DUP 1 # 130d               @ that below the x-axis line by line.
  113. PUT SUB PICT 3                    @
  114. ROLLD REPL -1                     @
  115.               STEP                @
  116.             ELSE 0                @
  117. 130                               @ If the x-axis is in the upper half  
  118.               FOR x               @ of the screen, generate the image   
  119. x R\->B lside 2 GET 2             @ from the x-axis to the bottom of    
  120. \->LIST PX\->C ITERATE            @ the screen.                         
  121.               NEXT                @ Start with the x-axis from left to  
  122. lside 2 GET 1 + B\->R             @ right.                              
  123. 63                                @                                     
  124.               FOR y               @                                     
  125. # 130d y R\->B 2                  @ Down left side to bottom of screen.
  126. \->LIST PX\->C ITERATE            @                                     
  127.               NEXT                @                                     
  128. 129 0                             @                                     
  129.               FOR x               @                                     
  130. x R\->B # 63d 2 \->LIST           @ Along bottom of screen right to     
  131. PX\->C ITERATE -1                 @ left.                               
  132.               STEP                @                                     
  133. 62 lside 2 GET 1 +                @
  134. B\->R                             @                     
  135.               FOR y               @                     
  136. # 0d y R\->B 2 \->LIST            @ Up right side to x-axis.
  137. PX\->C ITERATE -1                 @                     
  138.               STEP                @                     
  139. lside { # 130d                    @ Pass upper left and lower right     
  140. # 63d } MARIANI 0                 @ corners of rectangle to MARIANI.    
  141. lside 2 GET B\->R                 @                                     
  142.               FOR y               @ Map area below the x-axis to that   
  143. # 0d y R\->B 2 \->LIST            @ above the x-axis line by line.      
  144. PICT OVER PX\->C CONJ             @                                     
  145. C\->PX DUP 1 # 130d               @                                     
  146. PUT SUB PICT 3                    @                                     
  147. ROLLD REPL                        @                                     
  148.               NEXT                @                                     
  149.             END                   @                                     
  150.           \>>                     @                                     
  151.         ELSE 1 130                @                                     
  152.           FOR x x                 @ If the x-axis is not present, just  
  153. R\->B # 0d 2 \->LIST              @ ITERATE around the edges of the     
  154. PX\->C ITERATE                    @ screen to start.                    
  155.           NEXT 1 63               @ Along top.                          
  156.           FOR y                   @                                     
  157. # 130d y R\->B 2                  @                                     
  158. \->LIST PX\->C ITERATE            @ Down right side                     
  159.           NEXT 129                @                                     
  160. 0                                 @                                     
  161.           FOR x x                 @                                     
  162. R\->B # 63d 2 \->LIST             @                                     
  163. PX\->C ITERATE -1                 @ Along bottom                        
  164.           STEP 62 1               @                                     
  165.           FOR y                   @                                     
  166. # 0d y R\->B 2 \->LIST            @                                     
  167. PX\->C ITERATE -1                 @ Up left side                        
  168.           STEP {                  @                                     
  169. # 0d # 0d } {                     @ Pass upper left and lower right     
  170. # 130d # 63d }                    @ coordinates of the screen to        
  171. MARIANI                           @ MARIANI                             
  172.         END TICKS t               @                                     
  173. - B\->R 29491200 / f              @ Done! Calculate execution time.     
  174. STOF                              @ Restore flags.                      
  175.       \>>                         @                                     
  176.     \>>                           @                                     
  177.   MARIANI                         @                                     
  178. @ MARIANI is a recursive routine.
  179. @ It takes a rectangle (defined by its upper left corner (ulc) and lower right
  180. @ corner (rlc), divides it in two and iterates along the dividing line.  It
  181. @ then checks the border of each new rectangle in turn to determine whether all
  182. @ points on the border have the same status (either on or off).  If so, fill
  183. @ the rectangle.  If not, MARIANI calls itself with the coordinates of the
  184. @ smaller rectangle.  This contues until a rectangle is filled or until the
  185. @ current rectangle is less than 6 pixels in area in which case the state of
  186. @ each pixel is determined individually and no further division occurs.
  187.     \<< \-> ulc lrc               @
  188.       \<< lrc 1 GET               @                                     
  189. ulc 1 GET - lrc 2                 @                                     
  190. GET ulc 2 GET - *                 @                                     
  191.         IF # 6d <                 @
  192.         THEN ulc 1                @ If area less than 6 pixels,         
  193. GET 1 + B\->R lrc 1               @ determine state of each pixel in    
  194. GET 1 - B\->R                     @ rectangle individually.             
  195.           FOR x ulc               @                                     
  196. 2 GET 1 + B\->R lrc 2             @                                     
  197. GET 1 - B\->R                     @                                     
  198.             FOR y x               @                                     
  199. R\->B y R\->B 2 \->LIST           @                                     
  200. PX\->C ITERATE                    @                                     
  201.             NEXT                  @                                     
  202.           NEXT                    @                                     
  203.         ELSE lrc 1                @ If area greater than 6, divide      
  204. GET ulc 1 GET - lrc               @ rectangle in two.                   
  205. 2 GET ulc 2 GET -                 @                                     
  206.           IF <                    @                                     
  207.           THEN lrc                @ If the rectangle is taller than it  
  208. 2 GET ulc 2 GET + 2               @ is wide, determine line which       
  209. / DUP ulc 2 ROT PUT               @ divides it into upper and lower     
  210. SWAP lrc 2 ROT PUT                @ halves.  Set flag 4 to indicate     
  211. 4 SF                              @ this.                               
  212.           ELSE lrc                @                                     
  213. 1 GET ulc 1 GET + 2               @ If the rectangle is wider than it   
  214. / DUP ulc 1 ROT PUT               @ is tall, determine the line which   
  215. SWAP lrc 1 ROT PUT                @ divides it into left and right      
  216.           END \->                 @ halves.                             
  217. lines linet                       @                                     
  218.           \<<                     @ lines: start coordinate of line
  219.             IF 4                  @ linet: end coordinate of line       
  220. FS?C                              @                                     
  221.             THEN                  @                                     
  222. lines 1 GET 1 + B\->R             @                                     
  223. linet 1 GET B\->R 1 -             @ If flag 4 set, draw the horizontal  
  224.               FOR x               @ dividing line.                      
  225. lines 1 x R\->B PUT               @                                     
  226. PX\->C ITERATE                    @                                     
  227.               NEXT                @                                     
  228.             ELSE                  @                                     
  229. lines 2 GET 1 + B\->R             @                                     
  230. linet 2 GET B\->R 1 -             @                                     
  231.               FOR y               @ Otherwise draw the vertical one.    
  232. lines 2 y R\->B PUT               @                                     
  233. PX\->C ITERATE                    @
  234.               NEXT                @
  235.             END ulc               @
  236. linet CHECKBORDER                 @ Check the border of the left (or    
  237.             IF 1                  @ upper) rectangle (defined by        
  238. FS? 2 FS? OR                      @ corners ulc and linet).             
  239.             THEN                  @ If either flag 1 or flag 2 is set,  
  240. PICT ulc 1 GET 1 +                @ then all pixels around the border   
  241. ulc 2 GET 1 + 2                   @ of the rectangle have the same      
  242. \->LIST linet 1 GET               @ status.                             
  243. ulc 1 GET - 1 -                   @                                     
  244. linet 2 GET ulc 2                 @
  245. GET - 1 - BLANK                   @
  246.               IF 2                @
  247. FS?                               @ Fill black if flag 2 is set.
  248.               THEN                @ Fill white if flag 1 is set.        
  249. NEG                               @                                     
  250.               END                 @                                     
  251. REPL 1 CF 2 CF                    @                                     
  252.             ELSE                  @                                     
  253. ulc linet MARIANI                 @ If border pixels are not all one    
  254.             END                   @ state, call MARIANI with corners    
  255. lines lrc                         @ ulc and linet.                      
  256. CHECKBORDER                       @                                     
  257.             IF 1                  @ Check border of right (or lower)    
  258. FS? 2 FS? OR                      @ rectangle (defined by lines and     
  259.             THEN                  @ lrc).                               
  260. PICT lines 1 GET 1                @                                     
  261. + lines 2 GET 1 + 2               @                                     
  262. \->LIST lrc 1 GET                 @                                     
  263. lines 1 GET - 1 -                 @ Fill according to flag status as    
  264. lrc 2 GET lines 2                 @ above.                              
  265. GET - 1 - BLANK                   @                                     
  266.               IF 2                @                                     
  267. FS?                               @                                     
  268.               THEN                @                                     
  269. NEG                               @                                     
  270.               END                 @                                     
  271. REPL 1 CF 2 CF                    @                                     
  272.             ELSE                  @                                     
  273. lines lrc MARIANI                 @ Or call MARIANI with corners lines  
  274.             END                   @ and lrc.                            
  275.           \>>                     @                                     
  276.         END                       @                                     
  277.       \>>                         @                                     
  278.     \>>                           @                                     
  279.   CHECKBORDER                     @                                     
  280.     \<< \-> ulc lrc               @ Determine status of rectangle       
  281.       \<< ulc PIX? \->            @ border.                             
  282. state                             @ First determine state of ulc pixel.
  283.         \<< ulc 1 GET             @                                     
  284. 1 + B\->R lrc 1 GET               @                                     
  285. B\->R                             @                                     
  286.           FOR x ulc               @                                     
  287. 1 x R\->B PUT PIX?                @ Now check each pixel around border  
  288.             IF                    @ starting with top line and compare
  289. state \=/                         @ to state of ulc.  If they are the
  290.             THEN 5                @ same continue around border         
  291. SF lrc 1 GET B\->R                @ otherwise break out of the loop,    
  292. 'x' STO                           @ there is no need to check further.
  293.             END                   @ Set flag 5 to stop further          
  294.           NEXT                    @ checking.                           
  295.           IF 5 FC?                @                                     
  296.           THEN ulc                @
  297. 2 GET 1 + B\->R lrc 2             @ If flag 5 still clear check down    
  298. GET B\->R                         @ right side.                         
  299.             FOR y                 @                                     
  300. lrc 2 y R\->B PUT                 @                                     
  301. PIX?                              @                                     
  302.               IF                  @                                     
  303. state \=/                         @                                     
  304.               THEN                @                                     
  305. 5 SF lrc 2 GET B\->R              @                                     
  306. 'y' STO                           @
  307.               END                 @                                     
  308.             NEXT                  @                                     
  309.           END                     @                                     
  310.           IF 5 FC?                @                                     
  311.           THEN ulc                @ If flag 5 still clear, check left   
  312. 2 GET 1 + B\->R lrc 2             @ side.                               
  313. GET B\->R                         @                             
  314.             FOR y                 @                             
  315. ulc 2 y R\->B PUT                 @                             
  316. PIX?                              @                             
  317.               IF                  @                             
  318. state \=/                         @                             
  319.               THEN                @                             
  320. 5 SF lrc 2 GET B\->R              @                             
  321. 'y' STO                           @
  322.               END                 @                             
  323.             NEXT                  @                             
  324.           END                     @                                     
  325.           IF 5 FC?                @                                     
  326.           THEN ulc                @                                     
  327. 1 GET 1 + B\->R lrc 1             @ If flag 5 still clear, check        
  328. GET 1 - B\->R                     @ bottom.                             
  329.             FOR x                 @
  330. lrc 1 x R\->B PUT                 @
  331. PIX?                              @
  332.               IF                  @
  333. state \=/                         @
  334.               THEN                @
  335. 5 SF lrc 1 GET B\->R              @
  336. 'x' STO                           @
  337.               END                 @
  338.             NEXT                  @
  339.           END                     @
  340.           IF 5 FC?C               @
  341.           THEN                    @ If flag 5 clear then border all one
  342. state 1 + SF                      @ state.  Set flag 1 for all pixels   
  343.           END                     @ off.  Set flag 2 for all pixels     
  344.         \>>                       @ on.                                 
  345.       \>>                         @
  346.     \>>                           @
  347.   ITERATE                         @ ITERATE determines whether a        
  348.     \<< DUP DUP PIXON             @ coordinate is in the Mandelbrot set
  349. PIXOFF C\->R                      @ or not.  Toggle the current pixel   
  350. 0 \-> x y i                       @ as a visual aid in following        
  351.       \<< 0 0 0 0                 @ progress of program.                
  352.         WHILE SQ                  @ Initialize Z(0) as (0,0) and        
  353. SWAP SQ SWAP DUP2 +               @ duplicate.  Real and imaginary      
  354. 4 < 'i' INCR nITTR                @ parts are treated separately on the
  355. < AND                             @ stack.  The loop calculates         
  356.         REPEAT - x                @ Z(n+1)=SQ(Z(n))+c where c=x+iy.     
  357. + 3 ROLLD * DUP + y               @ It repeats until SQ(Z)>4 or nITTR   
  358. + DUP2                            @ iterations are completed.           
  359.         END 4 DROPN               @                                     
  360. x y R\->C i                       @                                     
  361.         IF nITTR \>=              @                                     
  362.         THEN PIXON                @ If nITTR iterations are completed,  
  363.         ELSE                      @ the coordinate is in the Mandelbrot
  364.           IF i 20 <               @ set.  Turn the pixel on.            
  365. 3 FS? AND                         @                                     
  366.           THEN i 2                @ If dithering is on, check to see    
  367. MOD                               @ whether the point diverged in less  
  368.             IF                    @ than 20 iterations.  If so and it   
  369.             THEN                  @ is odd, turn the pixel on.          
  370. PIXON                             @
  371.             ELSE                  @
  372. PIXOFF                            @
  373.             END                   @
  374.           ELSE                    @
  375. PIXOFF                            @
  376.           END                     @
  377.         END                       @
  378.       \>>                         @
  379.     \>>                           @
  380. END                               @
  381.